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1 . 


INTRODUCTION 


1.1 BACKGROUND AND INVESTIGATION OBJECTIVES 

The satellite microwave altimeters aboard GEOS-3 and 
SEASAT provided nearly global data on sea surface height (SSH) 
with respect to the reference ellipsoid. The altimeter data 
were used for a variety of oceanographic investigations, e.g. , 
Gulf-Stream dynamics (Ref. 1), estimation of mean sea surfaces 
(Refs. 2 and 3), and SSH mesoscale variability (Refs. 4 and 20). 

The feasibility of detecting and locating dynamic 
ocean currents using SSH data has been demonstrated with GEOS 
and SEASAT altimetry (e.g., Refs. 5 and 6). Proposed future 
satellite altimeter missions, such as the National Oceanic 
Satellite System (NOS3), would provide global SSH data for 
both oceanographic research and the generation of operational 
products (e.g., maps of geostrophic boundary currents, ring 
currents, and geostrophic velocity estimates). 

The purpose of this investigation is twofold: 

• Quantify the oceanographic information 
content of SEASAT altimeter SSH data 

• Develop and verify algorithms, suitable 
for use by NOSS , for automatically de- 
tecting and locating geostrophic ocean 
currents, eddy boundaries and rings, and 
for estimating geostrophic velocities 
from single tracks of satellite altimeter 
SSH data. 
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1.2 


TECHNICAL APPROACH 


The technical approach used in this investigation is 
based on three complementary techniques; 


• Autoregressive statistical modeling of 
the altimeter data 

• Matched- filter signal detectors for op- 
timally discriminating between noise and 
ocean-current signatures in the altimeter 
data 

• Kalman smoothing techniques for estimating 
oceanographic signal waveforms in the 
altimeter data. 


The application of these techniques is outlined in the follow- 
ing three sections. 

1.2.1 Study of the Oceanographic Inforroation Content 
of SEASAT Altimeter Data 


SEASAT altimeter data along nearly repeating ground 
tracks (repeat- track data) are used to study the oceanographic 
information content of the altimeter data. Geoid profiles are 
subtracted from each track of altimeter data. The statistical 
properties of the resulting residual time series are quantified 
by using autoregressive (AR) modeling to estimate the power 
spectrum of the residual data process. In addition, quantile- 
quantile plots are computed to characterize the probability 
distribution of the residuals. 

To reduce the geoid modeling errors to negligible 
magnitudes and to provide direct observations of the temporal 
variability of the oceanographic signals in the residual data, 
the data from repeating tracks are averaged, and the devia- 
tions of each track of data from this average are computed. 


1-2 



The statistics of the resulting difference data are analyzed 
by autoregressive modeling, and Kalman- smoothing techniques 
are used to estimate the waveforms in the data caused by meso- 
scale time-varying signal components. 

1.2.2 Development of an Ocean- Curt >;!.nt Detection Algorithm 

Matched- filtering techniques are used to develop an 
ocean-current detection algorithm that incorporates both de- 
terministic and statistical information about the ocean cur- 
rents, the altimeter data, and the estimated geoid profiles. 
This approach is flexible and suited for both boundary cur- 
rents and eddies. The matched- filtering approach is optimal 
with respect to reascoable models for the noise signals and 
the ocean-current signatures in the altimeter data. 

1.2.3 Development of a Geostrophic-Velocity 
Estimation A1 gor i t hm 


A bank of matched filters is used to detect the bound- 
ary current signature in the altimeter data and to compute 
maximum- likelihood estimates of the coefficients of a parametric 
model signature. The sea- surface slope is estimated from the 
parametric signature model, and the geostrophic equation yields 
an estimate for the cross- track component of the geostrophic 
current velocity. 


1.3 ORGANIZATION OF REPORT 

The report is organized as follows. The study of the 
oceanographic information content of SEASAT altimeter data is 
described in Chapter >. These results provide the foundation 
for the ocean-current detection algorithm, the development and 


1-3 



verification of which is presented in Chapter 3. Chapter A 
summarizes the development of the geostrophic-velocity estima- 
tion algorithm and its verification with SEASAT altimeter data 
in the western North Atlantic. The report concludes with a 
summary, conclusions, and suggestions for further study in 
Chapter 5. 



2 . 


OCEANOGRAPHIC INFORMATION CONTENT 


2 . 1 INTRODUCTION 

The approach for determining the oceanographic infor- 
mation content of SEASAT altimeter data starts with tne selec- 
tion of tracks for analysis from areas for which TASC has 
precise local geoid models. Geoid profiles are subtracted 
from each track of altimeter data to yield residual data sets, 
which consist of oceanographic signals as well as error signals 
caused by such factors as residual geoid modeling errors, orbit 
uncertainties, instrumentation noise, and errors in corrections 
for tides, barotropic pressure, and atmospheric influences. 

To study the signal components caused by time -varying 
oceanography, the residual data are analyzed along nearly re- 
peating tracks at 3-day intervals. Repeat-track analysis con- 
sists of first computing the point-wise ensemble mean of all 
the tracks in a given set. This ensemble mean is a time series 
that consists primarily of geoid modeling error together with 
the mean dynamic sea surface topography. The second step is 
to subtract this mean time series from each track in the set. 

The resulting difference data show how each track differs from 
the ensemble mean. These tracks of difference data give direct 
observations of the time-varying components of oceanography in 
the original altiiiOter data. 

The time-varying oceanographic signals in the differ- 
ence data are caused by both mesoscale and microscale variations 
in sea surface height. To determine the distribution of power 
with wavelength, the difference data are spectrally analyzed 
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by using autoregressive (AR) modeling techniques. The result- 
ing AR models yield an estimated power spectrum for each track. 
The differences among these spectra quantify the changes in 
the distribution of power in the oceanographic signals from 
one track to the next. The differences among these spectra 
may also be caused, in part, by time-varying residual error 
signals associated with the data corrections for atmospheric 
effects and instrument errors (as provided in the Geophysical 
Data Records). 

The power spectra provide estimates of the average 
power distribution in the difference data for wavelengths 
ranging from approximately 1500 km to 14 l<m. Based on this 
information, optimal Kalman smoothers are designed to extract 
the mesoscale (correlated) signal waveforms in each track of 
data. The optimal smoothers also compute theoretical rms 
values for the errors in the estimated mesoscale waveforms. 
These data provide direct quantitative measures th oceano- 
graphic information content of the SEASAT alti’»' data. 


2.2 SELECTED SEASAT ALTIMETER DATA 

Twenty-nine tracks of SEASAT altimeter data (1 sample/ 
second) were analyzed. These tracks were selected in the Gulf 
of Alaska and the North-Atlantic region south of Iceland, areas 
for which TASC has precise local 5-min gravimetric geoid esti- 
mates. Figures 2.2-1 and 2.2-2 show maps of these regions 
with 500-meter bathymetry contours. The selected SEASAT tracks 
are labeled with their REV numbers. 

Geoid profiles were computed along each track by using 
bilinear interpolation on precise 5-min local gravimetric geoid 
estimates. An example of a geoid profile is shown in Fig. 2.2-3 
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Figure 2.2-1 SEASAT Tracks and 500-Meter Bathymetry 

Contours for Gulf of Alaska 
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Figure 2.2-2 SEASAT Tracks and 500-Meter Bathymetry 

Contours for North Atlantic 
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Figure 2.2-3 Geoid Profile for SEASAT Rev. 1159D 

in North Atlantic 


for Rev. 1159D in the North Atlantic. The corresponding track 
of SEASAT altimeter data is shown in Fig. 2.2-4. Subtracting 
the geoid profile from the altimeter data produces the residual 
data shown in Fig. 2.2-5. The abrupt 0.4-meter feature located 
at 1050 km in Fig. 2.2-5 is caused by a seamount, while the 
offset and linear trend in this figure are attributed primarily 
to radial orbit uncertainty and long-wavelength geoid estimation 
error. The process of subtracting geoid profiles was carried 
out for all of the Revs, shown in Figs. 2.2-1 and 2.2-2. The 
resulting collection of SEASAT-minus-geoid residual data forms 
the basis for the analyses described in the following sections 
of this chapter. 
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Figure 2.2-5 SEASAT-Minus-Geoid Residual Data 

for Rev. 1159D in North Atlantic 
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2.3 


REPEAT -TRACK DATA 
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In the Gulf-of-Alaska region, nine of the satellite 
tracks (spanning 15 September - 7 October) were nearly repeating 
with an rms separation of 1.9 km. Eight more repeat tracks 
(spanning 16 September - 7 October) were available in the North- 
Atlantic region with an rms separation of 3.1 km. To analyze 
the time-varying oceanographic signals in these data, the long- 
wavelength errors caused by. radial orbit uncertainties were 
attenuated by subtracting a least-squares linear trend from 
each track of data. The resulting sets of residual data are 
shown in Figs. 2.3-1 and 2.3-2. For clarity, the plots are 
vertically offset from each other by 0.5 m. There are common 
features shared by all tracks in each set; these features are 
attributed to geoid estimation error and the mean sea-surface 
height signature. 

hiORTH-ATLANTIC REPEAT TRACKS 


R.730S* 



Figure 2.3-1 Ensemble of Detrended SEASAT-Minus-Geoid 

Residuals for Repeat Tracks in North Atlantic 
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Figure 2.3-2 Ensemble of Detrended SEASAT-Minus-Geoid 

Residuals for Repeat Tracks in Gulf 
of Alaska 


2.4 AUTOREGRESSIVE DATA ANALYSIS 

To determine the distribution of average signal power 
with wavelength (frequency), selected tracks of data were sta- 
tistically analyzed by using an autoregressive modeling tech- 
nique. This method of spectrum estimation is better suited 
for the analysis of individual tracks of data than the clas- 
sical method based on periodgrams. 

An autoregressive (AR) model of order p for a time 
series D(t), t = 1,..., N, is the difference equation 

P 

D(t) = Cj^D(t-k) + w(t), t = p+1,..., N (2.4-1) 
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where 




w(t) = residual noise 


The coefficients Cj^, k = 1,..., p, are chosen to minimize the 
mean-square residuals (COVAR algorithm, Ref. 9) 


\>’AR 


N 

, 1 ^ 
" N-p 


t=p+l 


w2(t> 


(2.4-2) 


If the AR model is appropriate for the process gener- 
ating the data, then the residuals w(t), t = p+1,..., N, are a 
sample of approximately white noise. It follows that the power 
spectral density (power spectrum) of the discrete-time process 
generating the data D(t) can be estimated as 


S„(F) = 


VAR 


- E 

k=l 


^k ® 


i2nFk 


(2.4-3) 


wh«re 


r?’ = dimensionless frequency (cycles/sample) (2.4-4) 


F = 0.5 is the folding frequency. The variance of the random 
process having the power spectrum S^(F) is 


variance 



S^(F)dF 


(2.4-5) 


A natural estimate for the power spectrum of the under- 
lying continuous-time process (of which the data D(t) are 
uniformly-spaced sample values) is 


V 


# 
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(2.4-6) 


S(f) = S^(f/fg)/fg. 



where 


f = data sampling frequency = 1.0 Hz 
s 

f = spectrum frequency in hertz 

The altimeter data D(t) are measured in meters, and it follows 

2 

that S(f) has the units m /Hz. The wavelength \ corresponding 
to frequency f is 


\ = v/f 


(2.4-7) 


where 


V - altimeter subtrack velocity = 6700 m/s (2.4-8) 

For each track of data, the best choice of the AR 
model order p is estimated by computing the Akaike information 
criterian (AIC) for a family of AR models, p = 0,..., N/20 
(Refs. 10-12): 

AIC = N log^(VAR) + 2p (2.4-9) 

The model for which AIC is smallest is chosen as the best AR 
model for the underlying process that generated the observed 
data D( t) . 

Several applications of this AR-modeling technique to 
the analysis of SEASAT data are described in the following 
sections of this chapter. 
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2.5 


PROBABILITY DENSITY FUNCTION OF RESIDUAL ALTIMETRY 


Selected tracks of SEASAT-minus-geoid residual altim- 
eter data were analyzed to detennine if their statistical be- 
havior is nearly Gaussian. Information about the probability 
distribution of the residual altimeter data is needed for de- 
signing optimal algorithms for detecting and estimating ocean 
currents and for computing the theoretical performance of these 
algorithms . 

To investigate the probability distribution of the 
altimeter data, the following approach was used: 


• Each track of residual altimetry was 
approximately whitened by applying a 
suitable linear transformation to the 
data. This linear transformation was 
determined from the minimum-AIC AR data 
model and was unique for each track of 
data. Mathematically stated, the data 
D(t) were transformed to the nearly white 
sequence w(t) by using the AR equation 

w(t) = D(t) - C^D(t-l) - ... - CpD(t-p) 

(2.5-1) 


where the coefficients Cj^, k = l,...,p, 

and the order p were determined as de- 
scribed in Section 2.4. 

• The whitened data w(t) were then analyzed 
to determine how nearly their empirical 
distribution function matched the normal 
(Gaussian) distribution. The use of 
whitened data is necessary because the 
statistical tests for normality are based 
on the assumption that the data are sta- 
tistically independent samples. 


The whitened data w(t) were listed in order from smal- 
lest to largest. These ordered values were then plotted versus 
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their theoretical deviation from the median, expressed in stan- 
dard deviations (quantiles) of the normal probability distribu- 
tion. The resulting graph is called a normal QQ plot, and (as 
depicted in Fig. 2.5-1) it shows departures of the data values 
from a normal distribution. Nearly Gaussian data yield a near- 
ly straight line, while short- or long-tailed empirical dis- 
tributions produce marked deviations from the straight line. 


R.7334S 



Figure 2.5-1 Example QQ Plots 


A typical QQ plot for residual SEASAT altimeter data 
is shown in Fig. 2.5-2; the plot indicates nearly Gaussian 
behavior. In contrast, Fig. 2,5-3 shows the worst-case (least 
Gaussian) result for all the tracks analyzed. The departures 
from Gaussian behavior are caused by atypical events such as 
the occurrence of data errors, extreme weather, and seamounts. 
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ORDERED DATA 















The conclusion reached from this analysis is that 
both Atlantic and Pacific residual altimeter data are nearly 
Gaussian in their statistical behavior. No systematic depar- 
ture from Gaussian behavior was detected. 


2.6 ANALYSIS OF TEMPORAL VARIABILITY IN REPEAT-TRACK DATA 

To study the time-varying signals in the repeat- 
track residuals, an ensemble-mean time series X was computed 
for each set of repeat- track residual data; 

m 

X(t) = i t=l,...,N (2.6-1) 

j=l 


where 

x.(t) = t-th datum in j-th track of residual data 
m = number of repeat tracks in set 
N = number of data in each track. 

The ensemble-mean series was then subtracted from 
each track of residual data to produce the difference-data 
time series D.. 

J • 


Dj(t) = x.(t) - X(t), 

j 1,..., m, t " 1,..., N 


( 2 . 6 - 2 ) 


The D.(t) time series (difference data) are direct observe- 

J 

tions of the time-varying oceanography in the original altimeter 
data. Because of the differencing operations, these difference 
data contain a negligible signal component caused by geoid 
estimation error. (The residual geoid estimation error in the 
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difference data has an estimated rms value of less than 0.3 cm 
based on the known rms separations of the repeat tracks and 
the known rms value of the first differences (X(t-l) - X(t)) 
of the ensemble-mean time series.) 

An example of SEASAT difference data is shown in 
Fig. 2.6-1 for Rev. 1159D in the North Atlantic. Because these 
data have a significant uncorrelated signal component, each 
data sample is plotted for clarity as a vertical spike. The 
rms value of these data is 6. A cm. The average rms value for 
the ensemble of eight tracks of difference data in the North- 
Atlantic region is 6.3 cm, the smallest rms value is A. A cm, 
the largest is 9.2 cm. In the Gulf-of-Alaska region the aver- 
age rms value for the nine repeat tracks of difference data is 
A. 2 cm, the smallest rms value is 3.0 cm, the largest is 6.5 cm, 
and the standard deviation of the rms values is 1.2 cm< 
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Figure 2.6-1 Difference Data for Rev. 1159D 

in North Atlantic 
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An AR spectrum estimate from the difference data along 
Rev. 1159D in the North-Atlantie region is shown in Fig. 2.6-2. 
This spectrum is representative of those computed for the other 
repeat tracks. The spectrum attains its largest values on the 
plateau at low frequencies (wavelengths greater than 1000 km) . 
There is a downward transition to another plateau at high fre- 
quencies (wavelengths shorter than 100 km). This high-frequency 
plateau is sometimes obscured by ripples in the spectrum esti- 
mate, as seen, e.g., in Fig. 2.6-2. These ripples may be caused 
in part by the inherent sampling variability of the spectrum 
estimator. The ripples may also be caused in part by high- 
frequency structure in the altimeter data. In either case, 
the overall structure of each AR spectrum estimate is repre- 
sented by a white-noise floor at high frequencies plus a 1st- 
order Markov spectrum that models the transition to the raised 
plateau at low frequencies. 
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Figure 2.6-2 Power Spectra for Difference Data Along 

Rev. 1159D in North Atlantic 
2-15 



Such a model spectrum is shown as the smoothly vary- 
ing sol-f * curve in Fig. 2.6-2, labeled "MWN MODEL SPECTRUM". 

The asymptotes of this model spectrum are shown as dashed 
lines. The lower horizontal asymptote is the spectrum of the 
white-noise component alone. The sloping asymptote is the 
high-frequency trend of the Markov component in the absence of 
white noise. The sum of the Markov noise and the white noise 
has the smooth spectrum shape shown as the MWN (Markov plus 
White Noise) model in Fig. 2.6-2. 

_2 

The corner frequency of 1.4 x 10 Hz marks the inter- 
section of the low-frequency plateau and the sloping asymptote 
of the Markov spectrum. The wavelength of this corner frequency 
is 480 km. 

With some tracks of data from the Gulf of Alaska, the 
low- frequency raised plateau is missing from the power spectrum 
because the Markov component is statistically insignificant; 
this situation is modeled by setting the variance of the Markov 
component equal to zero. 

The white-noise component of the difference data is 
attributed to uncorrelated altimeter instrument noise plus 
microscale oceanography. The Ist-order Markov component, in 
contrast, is caused by a correlated signal that is attributed 
primarily to mesoscale time-varying oceanography. The Markov- 
plus-white-noise (MWN) model for the difference data is impor- 
tant because it provides quantitative information about the 
oceanographic information content of the altimeter data. 

Data Model - The power spectrum of the MWN model has 
the following representation: 
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= standard deviation of Markov component (m) 
= corner frequency (hz) 
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^ww ” spectral density of white noise (m /Hz) 


The autocorrelation function of the Markov process is the in- 
verse Fourier transform of the Markov spectral density 
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where 

t = time shift parameter (s) 


The correlation time t^ of the Markov signal is defined as 

that value of the time shift for which R (t^) = R„. (0)/e: 

mm o mm ' 



( 2 . 6 - 6 ) 


The correlation distance d^ is computed by using the altimeter 
subtrack velocity v; 
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The variance of the white noise is computed for frequencies 
bounded by the folding frequency (half the data sampling 
frequency) : 


o 



f^/2 

-^s/2 


^ww^s 


( 2 . 6 - 8 ) 


where 


fg = 1.0 Hz for 1- sample/second SEASAT data.. 


For Rev. 1159D in Fig. 2.6-2, the rms value of the 
Markov component is = 4.3 cm, the correlation distance is 
d^ = 77 km, and the rms value of the white noise is o^ = 4.6 cm. 
The model parameters for all tracks of difference data are 


listed in Table 2.6-1 for the North Atlantic and in Table 2.6-2 


for the Gulf of Alaska. The average rms value of the Markov 
component in the North Atlantic is 4.0 cm as compared to 1.8 cm 
in the Gulf of Alaska. There is less difference between the 
average rms values of the white noise components; 4.8 cm in 
the North Atlantic as compared to 3.8 in the Gulf of Alaska. 
These rms levels are consistent with a recent study of meso- 
scale oceanographic variability described in Ref. 4. 


All of the analyzed repeat- track difference data were 
consistent with MWN models. As the data in Tables 2.6-1 and 
2.6-2 indicate, the parameter values of the MWN models vary 
significantly from track to track in both the North-Atlantic 
and Gulf-of-Alaska regions. The variability of model param- 
eters is depicted in Fig. 2.6-3, which shows the average of 
the spectrum aymptotes for the North-Atlantic data together 
with their minimum, maximum, and rms variations. 
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TABLE 2.6-1 

MODEL PARAMETERS FOR REPEAT-TRACK DIFFERENCE DATA 

NORTH ATLANTIC 


REV. NO. 

DAY 

RMS OF MARKOV 
COMPONENT 
(cm) 

CORRELATION 

DISTANCE 

(km) 

RMS OF 
WHITE NOISE 
(cm) 

1159D 

0 

4.3 

77 

4.6 

1202D 

3 

3.5 

52 

3.1 

1245D 

6 

4.9 

137 

5.4 

1288D 

9 

5.5 

10 

7.4 

1331D 

12 

2.6 

15 

4.0 

1374D 

15 

2.9 

49 

3.3 

1417D 

18 

3.6 

28 

5.4 

1460D 

21 

4.4 

21 

5.1 

Average Value 


4.0 

49 

4.8 

Standard 

Deviation 


1.0 

42 

1.4 

Maximum Value 


5.5 

137 

7.4 

Minimum Value 


2.6 

10 

3.1 


The observed variability indicates that the time- 
varying signals in repeat-track data are best modeled as non- 
station, ary from track to track. Along single tracks the 
difference data are consistent with stationary models (for 
track lengths of 1300 km in the North Atlantic and 1000 km in 
the Gulf of Alaska). The parameter values of these stationary 
models sometimes change significantly over the three-day inter 
val between successive repeat tracks. 
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TABLE 2.6-2 

MODEL PARAMETERS FOR REPEAT-TRACK DIFFERENCE DATA 

GULF OF ALASKA 


REV. NO. 

DAY 

RMS OF MARKOV 
COMPONENT 
(cm) 

CORRELATION 

DISTANCE 

(km) 

RMS OF 
WHITE NOISE 
(cm) 

1155A 

0 

2.6 

21 

4.7 

1198A 

3 

3.3 

30 

3.9 

1241A 

6 

0 

NA 

4.0 

1284A 

9 

0 

NA 

3.0 

1327A 

12 

1.9 

69 

3.1 

1370A 

15 

2.4 

7 

2.2 

1413A 

18 

1 . 6 

44 

2.9 

1456a 

21 

0 

NA 

3.5 

1499A 

24 



4.1 

8 

5.1 

Average Value 

1.8 

30 

3.6 

Standard 

Deviation 


1.5 

24 

0.9 

Maximum Value 

4.1 

69 

5.1 

Minimum Value 

0 

7 

2.2 

NA = Not Applicable 




Kalman Smoothing - An important question bearing on 
the information content of the altimeter data is 


• How accurately can the signal waveforms 
of the mesoscale components in the dif- 
ference data be separated from the addi 
tive white noise? 
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Figure 2.6-3 Sunuaary of Model Power Spectra for 

Difference Data in North Atlantic 


To answer this question, a Kalman fixed- interval 
smoother (Ref. 13) was designed for each track of difference 
data in the North Atlantic. Each smoother was optimal for 
estimating the mesoscale (Markov) waveform in a particular 
track of altimeter data. The waveform estimates were optimal 
in the sense that the estimation error variances were mini- 
mized, given the data model derived from AR analyses for each 
track of data. 

The smoothed estimates of the mesoscale waveforms in 
the North-Atlantic data are shown in Fig. 2.6-4. Beside each 
waveform is the day number for the track and the ms error of 
the smoothed estimate as computed by the Kalman smoother for 
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that track. The theoretical rms errors of the estimates range 
from 1.3 cm to 2.5 cm, with an average of 1.9 cm. By comparing 
Fig. 2.6-1 (which shows the raw difference data for day 0> 
with the smoothed waveform at the top of Fig. 2.6-4, one sees 
how effectively the smoother suppresses the white noise in the 
data. 



• tM IMt 

DISTANCS 

Figure 2.6-4 Smoothed Estimates of Time-Varying 

Mesoscale Signals in Repeat-Track Data 
from North Atlantic 

The smoothed waveforms in Fig. 2.6-4 show statistically 
significant peak-to-peak differences from one track to the 
next. A good example of this variability is provided by the 
estimates for days 6, 9, and 12. The smoothed estimate for 
day 6 has a 19-cm peak-to-peak variation, while the day-9 and 
day-12 estimates are not significantly different from zero 
over most of their lengths. 
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2.7 TEMPORAL VARIABILITY OF GDR DATA CORRECTIONS 

The estimated mesoscale waveforms in Fig. 2.6-4 may 
be affected by errors associated with the data corrections 
provided in the Geophysical Data Record (GDR) for the altimeter 
instrument and atmospheric effects. This potential error source 
was investigated by analyzing the sum of the GDR instrument 
corrections plus the atmospheric corrections for the set of 
repeat tracks in the North Atlantic. 

The method of analysis followed the same pattern as 
the analysis of differenced altimeter data described earlier 
in this chapter. The analysis technique consisted of three 
steps : 

• The time series of GDR corrections for 
all the repeat tracks were averaged to 
form an ensemble-mean time series for 
the set of eight tracks 

• This ensemble-mean time series was sub- 
tracted from each track of correction 
data to form a time series of difference 
data for that track 

• The power spectrum was estimated for 
each track of difference data (one for 
each repeat track) by subtracting a least- 
squares linear trend from each time series 
and then using the autoregressive (AR) 
modeling technique described in Section 2.4. 

The difference data contain the time-varying components 
of the correction signals, and the power spectra measure the 
distribution of power in these signals as a function of fre- 
quency. The rms magnitudes of the time-varying correction 
signals varied from track to track: the minimum was 0.2 cm; 

the maximum was 1.6 cm; and the mean was 0.8 cm. The standard 
deviation about the mean was 0.4 cm. These rms magnitudes are 
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Gmail compared to the nns values found for the time-varying 
mesoscale signals in the sea surface height data. 

The estimated power spectra for the time-varying cor- 
rection signals were compared with the corresponding spectra 
for the time-varying signals in the sea-surface height (SSH) 
data. The comparison showed that the time-varying SSH data 
typically have more than ten times the power of the time- 
varying correction data for wavelengths shorter than 1000 km. 

If the assumption is made that the errors in the GDR 
corrections for Instrument and atmospheric effects are smaller 
than the corrections themselves, then the analysis of the cor- 
rection data supports the following conclusion: in the nearly 

repeating tracks of altimeter data used for this study, the 
instrvunent and atmospheric effects are significantly smaller 
than the observed temporal variability. 


2.8 CONCLUSIONS ON INFORMATION CONTENT 


Based on the analyses of North-Pacific and North- 
Atlantic residual SEASAT altimeter data, the following conclu 
sions are reached: 


• Residual altimeter data have the statis- 
tical properties of correlated Gaussian 
random noise 

• Residual data are statistically modeled 
as the sum of white (uncorrelated) noise 
plus a Markov (mesoscale) process 

• Mesoscale signal components in the resid- 
ual data are attributed primarily to 
geoid estimation errors and mean sea- 
surface height signatures 
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• Each track of repeat-track difference 
data has the statistical properties of a 
Gaussian process consisting of white 
noise plus a mesoscale lst>?order Markov 
process. The rms values of the white 
and mesoscale signal components often 
vary significantly over 3-day intervals 
within each set of repeat tracks. The 
average rms values of these signal compo- 
nents (and sample standard deviation o 
of the rms values) are listed for the 
two study regions: 

RMS WHITE-NOISE LEVELS 

North Atlantic: 4.8 cm (o = 1.4 cm) 

Gulf of Alaska: 3.6 cm (a = 0.9 cm) 

RMS MESOSCALE SIGNAL LEVELS 

North Atlantic: 4.0 cm (o = 1.0 cm) 

Gulf of Alaska: 1.8 cm (o = 1.5 cm) 

• Kalman smoothers for the difference data 
provide minimum-variance estimates of 
the time -varying mesoscale signal wave- 
forms. The nonstationary data require a 
unique smoother for each repeat track, 
and the rms estimation accuracy varies 
from track to track. Average rms estima- 
tion accuracy for North-Atlantic data is 
1.9 cm. 

• The mean- square temporal variability of 
instrumental and atmospheric corrections 
for the altimeter data (Sept. 1980 GDR 
tape) is less than 10 percent of the 
variability in the difference data for 
wavelengths < 1000 1cm. 


2.9 QUESTIONS FOR FURTHER STUDY 

The results of this study raise the following ques- 
tions for future study: 

• How much of the observed mesoscale vari- 
ability in the repeat-track difference 
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data is correlated with variations in 
the significant wave height on the GDR 
tape? 

• If there is significant correlation, how 
does it vary with wavelength? 


The answers to these questions may be found on a track-by- track 
basis via multi-channel AR modeling techniques for estimating 
spectral coherence. 


2-26 


3. 


OCEAN -CURRENT DETECTION ALGORITHM 


3 . 1 INTRODUCTION 


Approach - The approach to developing an automatic 
algorithm for detecting specific ocean-current signatures in 
single tracks of altimeter data is based on the results o£ the 
study of the oceanographic information content of SEASAT data 
in Chapter 2. The approach consists of three steps; 


• The residual altimeter data are modeled 
as samples of Gauss-Markov noise plus 
the possible occurrences of specific 
ocean-current signatures (boundary cur- 
rents, cold-core rings, etc) 

e The correlation structure of the noisy 

signal components in the data are identi- 
fied for each track of data with a data- 
adaptive AR modeling algorithm 

• The hypotheses that specific current 
signatures are present in the data are 
tested by using optimal matched filters 
and threshold detectors. 


Algorithm Structure - The ocean-current detection 
algorithm processes single tracks of residual satellite altim- 
eter data and yields the following outputs: 


• Detected locations of specified ocean- 
current signatures along the satellite 
subtrack 

• Estimated amplitudes of the detected 
signatures 
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• Estimated rms errors for the locations 
and amplitudes of detected signatures 

• Expected number of false alarms. 


The residual altimeter data are inputs to the algo- 
rithms and are computed from raw altimeter data in three steps 
by 


• Applying corrections for known error 
sources 

• Interpolating the data through intervals 
in which the data are in serious error 
(e.g. , outliers) 

• Subtracting an estimated gravimetric 
geoid profile along the satellite 

s$uui.i.aCiv.* 

The resulting residual data are noisy measurements of 
the dynamic sea- surface height. The characteristics of the 
noise in these data depend on the noise in the raw altimeter 
data and on the accuracies of the error corrections and the 
geoid profiles. The detection algorithm exploits both the 
statistical properties of the noise in the residual data and 
the known average properties of ocean-current signatures in 
the altimeter data. For specified models of the noise and 
oceanographic signature, the algorithm maximizes the probabil- 
ity of detection at a specified probability of false alarm and 
minimizes the rms errors of the estimated current signature 
parameters . 

As depicted in Fig. 3.1-1, the detection algorithm 
consists of four subalgorithms that perform separate functions. 

• AUTOREGRESSIVE MODELING - The residual 
data are analyzed to detennine a sto- 
chastic autoregressive (AR) model for 
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Figure 3.1-1 Structure of Data-Adaptive Current 

Detection Algoritlim 


the process that generated the data. 

The order of the AR model is selected to 
minimize the Akaike information cri- 
terion (AIC). 

• MATCHp-FILTER DESIGN - The AR model, 
together with a user-specified generic 
ocean-current signature, are used to 
compute the impulse response of the op- 
timal matched filter for detecting and 
locating the generic signature in the 
noisy residual data. 

• MATCHED FILTER - The impulse response of 
the matched filter is convolved with the 
residual altimeter data to compute a 
sequence of sufficient statistics for 
the threshold detector. 

• THRESHOLD DETECTOR - The detector compares 
the sufficient statistics with a threshold 
value that is chosen to yield a specified 
false-alarm rate. A detection occurs 
when the statistic exceeds the threshold. 
The estimated location of the detected 
signature is given by the location of 

the local maximum of the statistic. 


OCTECTION 

DATA 
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Overview - This chapter describes the dev<ulopment and 
verification of the detection algorithm. Section 3.2 describes 
the SEASAT data selected for the verification of the algorithm. 
Generic ocean-current signature models are presented in Sec- 
tion 3.3. The results of processing SEASAT data to detect the 
Gulf Stream and cold-core ring currents in the western North 
Atlantic are discussed in Sections 3.4 and 3.5. After discus- 
sing the verification results in Section 3.6, the mathematical 
basis of matched-filter detection is discussed in Sections 3.7 
and 3.8. Algorithm performance is discussed in Sections 3.9 
and 3.10. The chapter concludes with a summary in Section 3.11 
and a discussion of areas for future study in Section 3.12. 


3.2 DATA SELECTED FOR DEMONSTRATING CURRENT- 
DETECTION ALGORITHM 

Four tracks of SEASAT altimeter data from the western 
North Atlantic were selected to verify the performance of the 
detection algorithm. The tracks are Revs. 234A, 277A, 478A, 
and 572D. These tracks crossed cold ring No. 4 as depicted in 
Fig. 3.2-1. The motion of this ring current had been tracked 
•using a variety of oceanographic data (Ref. 6). These four 
tracks also intersected the Gulf Stream and therefore provided 
data to verify the algorithm performance with respect to the 
detection of both boundary currents’ and eddy currents. 

The tracks of data used to verify the algorithm ex- 
tend beyond the latitude limits shown in Fig. 3.2-1; the test 
data extend from 25 deg to 39 deg north latitude. 
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SEASAT TRACKS AND COLD RING N0.4 
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Figure 3.2-1 Selected Data for Demonstrating 

Current-Detection Algorithm 


3.3 GENERIC OCEAN-CURRENT SIGNATURE MODELS 


This section describes two parametric families of 
ocean-current signatures. The first family is used for design- 
ing matched filters to detect warm-core and cold-core current 
rings. The second family is intended for detecting boundary 
currents, such as the Gulf Stream, and for estimating geo- 
strophic current velocities. 


Ring-Current Signatures - A family of generic altim- 
etric signatures is described for modeling the dynamic sea- 
surface features caused by cold-core and warm-core current 
rings. The sea-surface height H(x) at radial position x with 
respect to the ring's center is modeled as 
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ORIGINAL PAGE IS 
OF POOR QUALITY 

H(x) = - D exp (-9.21 (x/W)^) 

D = signature depth 

W = signature width (10 percent depth) 

D is positive for cold rings and negative for warm rings. 

This parametric model has a simple mathematical form and ap- 
pears to be in reasonable agreement with available data on 
ring signatures (e.g. , Refs, 6-8). . 

The width W is defined as the diameter at which the 
signature is 10 percent of its central value: 

H(W/2) = H(0)/10 (3.3-2) 

Equation 3.3-1 has the form of a Gaussian probability 
density. Therefore, these are referred to as Gaussian ring 
signatures. An example of a Gaussian ring signature is shown 
in Fig. 3.3-1, where the central depth is 0.5 meter and the 
width is 150 km. 


GAUSSIAN COLD-RING MODEL n-7Mia« 



Figure 3.3-1 Gaussian Cold-Ring Signature, 

Depth = 0.5 m. Width = 150 km 
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The tangential current-velocity distribution implied 
by a ring signature may be computed by setting the radial slope 
of the sea surface equal to the sum of the horizontal Coriolis 
acceleration and the centrifugal acceleration divided by the 
acceleration of gravity: 

dH(x) _ f v(x) + v^(x)/x /o o_o\ 

dx “ g vo.o-o; 


f = 2fi sin(|) = Coriolis parameter 
n - earth's rotational velocity 
(]) = latitude 

v(x) = tangential current velocity 
g = acceleration of gravity. 


The geostrophic velocity component is 

Vg<x> = f (3.3-A) 

Solving Eq. 3.3-3 for the total current velocity v(x) yields 


v(x) 


X f 

T" 




AV (x) 
~xt 


- 1 


(3.3-5) 


For a Gaussian 0.5-m 150-km ring signature at 45-degrees 
latitude, Eqs. 3.3-4 and 3.3-5 yield the velocity distributions 
shown in Fig. 3.3-2. The geostrophic approximation is seen to 
over-estimate the maximum velocity by approximately 0.1 m/s 
(18 percent) . 


Boundary Current Signatures - A family of generic 
altimetric signatures for boundary currents is defined with 
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Figure 3.3-2 Tangential Current-Velocity Distributions in 

Gaussian Cold Ring. Crosses = geos trophic 
approximation; circles = geostrophic approx- 
imation with centrifugal correction. 


the aid of Fig. 3.3-3, which depicts a satellite subtrack cros 
sing a current at angle 0. At position x along the subtrack, 
the dynamic height H(x) is modeled with the hyperbolic tangent 
function: 


H(x) = -(A/2) tanh(3 x sin0/W^) (3.3-6) 

A = amplitude of dynamic height change 

0 = track angle with respect to current velocity 

W = width of current (90 percent height change) 
c 
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Figure 3.3-3 Geometry of Geos trophic Current 

and Satellite Subtrack 


As shovm in Fig. 3.3-4, the tanh model signature pro- 
vides a reasonable fit to the SEASAT altimeter data for Rev. 
234A where the subtrack crosses the Gulf Stream centered at 
sample number 233. The tanh model provides similar fits to 
other tracks of SEASAT data; further examples are discussed in 
Section 4.3. 


n«037s 


GULF-STREAM SIGNATURE 



Figure 3.3-4 Comparison of Tanh Model Signature 

with SEASAT Data Rev. 234A 
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The along- track slope of the tanh signature is 
dH(x) _ 3 A sin0 u2 

- • ~~r^~ 

For the coordinate system in Fig. 3.3-3, this slope is related 
to the cross-track component V^(x) of the geostrophic velocity 
as follows 


3 X sin8 
W. 


(3.3-7) 


V^(x) = - f (3.3-8) 

For tracks that cut across the current, the geostrophic 
velocity profile V (x) along the subtrack is proportional to 
the cross- track velocity V (x) 


V tv's = V /v^/c^nA 
'g'**' ’c'**'' 




For the signature slope given by Eq. 3.3-7, the geo- 
strophic velocity profile is therefore 


V <x) = g-g § sech^ 


3 X sind 


W, 


(3.3-10) 


The signature amplitude parameter A is proportional 

to the maximum geostrophic velocity V (0) 

g 

. N 

2 f W 

A = ^ g g Vg(0) (3.3-11) 


For tracks that intersect the current, the signature width W 

s 

is proportional to the current width 

= W /sine (3.3-12) 
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Typical model parameters for the Gulf Stream in the 
western North Atlantic are an amplitude of A = 1 m, a maximum 
geostrophic velocity of V„(0) = 2 m/s, and a latitude of 

Q 

(j) = 45 deg. From Eq. 3.3-11 the current'^^s width is = 71 km. 
For a nominal track crossing angle of 60 deg, the along- track 
width of the signature is W_ = 82 km. Figure 3.3-5 shows the 
dynamic sea-surface height signature for these parameter values, 
while Fig. 3.3-6 depicts the geostrophic velocity profile im- 
plied by the height signature. 


TANH GULF-STREAM MODEL 



Figure 3.3-5 Dynamic Sea-Surface Height Signature 





Figure 3.3-6 
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Geostrophic Velocity Profile 
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3.4 


DEMONSTRATIONS OF GULF-STREAM DETECTION 


This section presents four examples of using the algo- 
rithm (Section 3.1) to detect automatically the Gulf-Stream 
signature in SEASAT data from which Marsh-Chang geoid estimates 
(Ref. 14) were subtracted. The tanh model signature (Eq. 3.3-6 
with Wg = W^/sin0) was used with the signature width Wg ranging 
from 35 km to 150 km. 

Example 1 - As depicted in Fig. 3.4-1, the residual 
SEASAT data for Rev. 572D are used as input to the detection 
algorithm. The algorithm is optimized in this example for the 
tanh model signature having a width of 60 km. The scaled out- 
put of the matched- filter is shown with the detection threshold 
chosen for a false-alarm rate (FAR) of 0.01 false alarm per 
megameter of track length (A/Mm). 

The filter output is seen to rise above the threshold 
and reach a local maximum at data sample number 64. This is a 
correct detection of the Gulf Stream. There are no other cross- 
ings of the threshold, which means that there are no false 
alarms. The algorithm successfully discriminates against the 
signature in the altimeter data due to cold ring number 4. 

Figure 3.4-2 shows the same input data for Rev. 572D 
along with the resulting scaled outputs of three matched fil- 
ters. Each filter is optimized for a different signature width 
(40 l<m, 60 km, and 150 km). All filters correctly detect the 
Gulf Stream with no false alarms, which demonstrates that the 
algorithm performance is insensitive to the choice of model 
signatures. The detection threshold in each case is chosen to 
yield FAR = 0.01 A/Mm. The 60-km filter produces the largest 
Gulf-Stream response, and for that reason its estimate of the 
current's location (sample nvimber 64) is the most likely of 
the three estimates available. The theoretical basis for using 
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Figure 3.4-1 Demonstration of Gulf-Stream Detection Using 

SEASAT Data Rev. 572D 
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the filter that produces the largest scaled output is explained 
in Section 3.7. The autoregressive (AR) model for the data is 
computed by the detection algorithm operating on the last 
230 samples of the altimeter data; this avoids Gulf-Stream 
bias in the modeling of the background noise in the altimeter 
data. 


The output of each matched filter is scaled so that 
the random noise in the output signal has a standard deviation 
of unity. Therefore, the detection threshold of the 60-km 
filter, which is labeled 3.0, corresponds to an output level 
that is three sigmas (three standard deviations) away from the 
mean. The three filters differ in their dynamics, as indicated 
by the smoother output waveforms from the filters optimized 
for the wider signatures. Because of th^se differing dynamics, 
the fixed false-alarin rate leads to slightly different detec- 
tion thresholds (when expressed in standard deviations) for 
the three filters. 

Example 2 - Figure 3.4-3 depicts the results of proc- 
essing the residual SEASAT data for Rev. 478A with three matched 
filters optimized for widths of 35 km, 100 km, and 150 km; the 
detection thresholds are chosen to yield FAR = 0.01 A/Mm. In 
this example, the Gulf Stream is located at the right end of 
the data set because Rev. 478A is an ascending track in which 
the Gulf Stream occurs at the end of the time series. All 
three filters are successful in detecting the Gulf Stream, and 
there are no false alarms. The AR model for the data is gen- 
erated by the algorithm using the first 215 samples of data to 
avoid the Gulf -Stream signature. Experience shows that the 
algorithm performs well even if the Gulf-Stream signature is 
included in the data used to develop the AR model. The primary 
consequence of using all the data for AR modeling is that the 
detection thresholds are moved to slightly higher values, because 
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4-3 Demonstration of Gulf-Stream Detection 
SEASAT Rev. 47 8A 
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the algorithm then underestimates the ease with which the Gulf- 
Stream signature can be detected. 

Example 3 - Figure 3.4-4 is an example of Gulf-Stream 
detection using Rev. 277A. This Rev. is interesting because 
it contains erroneous data caused by the island of Bermuda. 

The outputs of two matched filters are shown, one matched to a 
narrow 35-km signature and the other matched to a 60-km signa- 
ture; these filters use an AR noise model based on the first 
225 data samples. Both filters correctly detect the Gulf Stream, 
and there are no false alarms. The data errors due to Bermuda 
cause a relatively large response from the narrow 35-km filter, 
which is the more sensitive filter to narrow waveforms in the 
data. This example illustrates that the 60-lon matched-filter 
detector is robust against naturally occurring localized data 
errors . 


Example 4 - This last example uses SEASAT data from 
Rev. 234A as shown in Fig. 3.4-5. The data set contains an 
isolated data error and the edge of a warm-ring signature. 

The scaled outputs of two matched filters are shown with detec- 
tion thresholds corresponding to FAR = 0.01 A/Mm. The filters 
are matched to signature widths of 35 km and 100 km, and the 
AR noise model is based on the first 225 data. Both filters 
correctly detect the Gulf Stream without any false alarms. 

The data error causes a relatively large output from the 35-km 
filter, but has little effect on the output of the 100-km fil- 
ter. Both filters discriminate well against the signature of 
warm ring S . 
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Figure 3.4-5 
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3.5 


DEMONSTRATIONS OF COLD-RING DETECTION 


This section presents four examples of using the ^ie- 
tection algorithm to detect automatically the signatures of 
cold ring mjimber 4 in SEASAT data from which Marah-Chang geoid 
estimates (Ref. 14) were subtracted. The data sets used here 
are the same as the ones used in Section 3.4. The Gaussian 
model signature (Eq. 3.3-1) is used in the algorithm with two 
values for the signature width: 150 km and 300 km. 


Example 1 - Figure 3.5-1 shows the SEASAT data for 
Rev. 572D and the outputs of two matched filters. The algo- 
rithm designed these filters for signature widths of 150 km 
and 300 km by using an AR noise model based on the last 230 
data (to avoid the Gulf-Stream signature). The detection 
thresholds correspond to FAR = 0.1 A/Mra. This false-alann 
rate is ten LXIueS iSrgei. uiiau uuc vaxuc xn 


for Gulf-Stream detection. The larger FAR is used in this 
example because cold-ring signatures are more difficult to 
detect than boundary- current signatures; the ring signatures 
are more difficult to discriminate against the ambient noise 


signals in the data. 


Cold ring number 4 is correctly detected by the 300-km 
filter. The 150-km filter misses the ring because the model 
signature is too narrow. Both filters produce false alarms 
due to the Gulf Stream, but these are the only false alarms. 
Experience shows that the Gulf Stream consistently causes false 
alarms when the algorithm is attemptitig to detect cold rings. 
The conclusion to be reached is that the data must be processed 
in two passes. The first pass is used to detect the possible 
occurrence of boundary-current signatures. If a boundary- 
current signature is detected, then the data set is truncated 
to exclude the detected current during the second pass. The 
second pass is then used to detect the possil le occurrence of 
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Figure 3.5-1 


MATCHED-FILTER OUTPUTS 



eddy currents. This strategy requires a small amount of com- 
putation (0.5 sec per pass using an IBM 4341 computer) and 
prevents false alarms due to boundary currents. This two-pass 
technique is used in the following examples. 

Example 2 - Figure 3.5-2 shows SEASAT data from 
Rev. 478 a, which are used as input for two filters matched to 
signatures having widths of 150 km and 300 km. The algorithm 
implemented these filters using an AR model based on the first 
140 data. This selective use of data leads to the most accu- 
rate detection thresholds to achieve a specified false-alarm 
rate. When the algorithm uses all the data to compute the AR 
noise model (as it likely would in operational data analysis), 
the detection thresholds are raised slightly because the algo- 
rithm then overestimates the likelihood of false alarms. 

The 300-km filter correctly detects cold ring number 
4 at a false-alarm rate of 0.03 A/Mm. In contrast, the 150-km 
filter misses the ring even with the larger FAR of 0.1 A/Mm, 
which demonstrates the ability of the algorithm to discriminate 
between narrower and wider current signatures. There are no 
false alarms. 

Example 3 - Figure 3.5-3 shows the SEASAT data from 
Rev. 277A, which contains data errors caused by the island of 
Bermuda and a signature due to cold ring number 4. The data 
are processed by two filters matched to cold-ring signatures 
having widths of 150 km and 300 km. The algorithm implemented 
the filters using an AR noise model based on the first 134 data 
to avoid any bias from the Bermuda signature. The 300-km filter 
correctly detects cold ring number 4 at a false-alarm rare of 
0.2 A/Mm. The depression in the altimeter data to the left of 
Bermuda also causes an alarm in the output. This may be a 
false alarm; no data are available for certifying the origin of 
the signature. The 150-km filter misses cold ring number 4. 
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Figure 3.5-3 
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Example 4 - Figure 3.5-4 shows the SEASAT data from 
Rev. 234A; the isolated data error has been corrected by sim- 
ple interpolation. The edge of cold ring number 4 is visible 
and spans approximately 200 km in the data. To detect this 
signature, the data are processed twice with the detection 
algorithm, once to implement a detection filter matched to a 
signature width of 150 km, and again for a filter matchr'.d to a 
width of 300 km. For this example, the algorithm uses all 
220 data in computing the AR noise model. 

The scaled outputs of the filters show that the cold 
ring is correctly detected by both filters at a false-alarm 
rate of 0.5 A/Mm. There are no false alarms with the 150-km 
filter. Whether the alarm in the middle of the 300-km filter 
output is false cannot be certified because surface truth data 
are not available. 


3.6 DISCUSSION OF DEMONSTRATION RESULTS 


The verifications of algorithm performance described 
in Sections 3.4 and 3.5 lead to several conclusions: 


• The required computing time is modest: 
less than 0.5 sec per track for a one- 
signature (one-filter) analysis using an 
IBM 4341 computer 

• The Gulf Stream is easily detected at a 

low false-alarm rate: a single model 

signature is adequate; and FAR < 0.01 A/Mm 

• Cold rings are detected with some risk 
of false alarms. Processing the altim- 
etry with two or three filters (to search 
for two or three model signatures) in- 
creases the likelihood of detection for 

a specified false-alarm rate. The detec- 
tion of cold ring number 4 typically 
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occurs at false-alarm rates in the range 
of 0.05 to 0.5 A/Mm 

• The detection results are insensitive 

to localized data errors, which suggests 
that special preprocessing to correct 
erroneous data is unnecessary 

• A two-pass technique for detecting bound- 

ary currents and ring currents is recom- 
mended for processing individual tracks 
of residual altimeter data: the first 

pass is used to detect possible boundary 
currents with a single model signature 
(one filter). If a boundary current is 
detected, the data set is truncated to 
remove the detected current signature. 

In the second pass , the truncated data 
set is processed with two or three fil- 
ters, each matched to a different-width 
ring signature (e.g., 300, 200, 100 km). 
At each position along the track, the 
threshold detector uses the filter with 
the largest scaled output value. 

• Since the model signatures for warm ring 
currents are the negatives of correspond- 
ing cold-ring signatures, a single filter 
serves to detect both cold and warm rings 
when a two-sided threshold detector is 
used. For example, if the filter is 
matched to a cold-ring signature and the 
detector uses a threshold TH > 0 , then 
corresponding warm rings can be mapped 

by using a detector that flags those 
excursions of the filter output that 
fall below the negative threshold 
value -TH. 


3.7 MATCHED- FILTER DETECTION THEORY 

The theory of optimal matched filters for detecting 
deterministic signatures in additive colored noise is discus 
sed in several text books (e.g., Refs. 15, 16, 18). The key 
results of that theory for ocean-current detection are summa 
rized in the following. 
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Statement of the Detection Problem - The problem of 
detecting ocean-current signatures in residual altimeter data 
is formalized as follows. 

GIVEN! D(t) = time series of residual altimeter data 

m(t> = ocean-current signature time series 

N(t) = stationary Gaussian noise model for 

residual altimeter data that are free 
of current signature m(t) 

T = specified time (location) in the data D(t) 

Ag = unknown signature amplitude scale factor 

Ep = hypothesis that ID(t) = N(t) + A_m(t-T) 

^ with Ag 7 ? 0 ® 

Hq = null hypothesis that D(t) = N(t) 

FIND: An optimal decision rule for correctly choosing 

between hypotheses Hq and Ej,; and an optimal 

estimate of the amplitude Ag when Ej, is chosen. 

OPTIMALITY: Maximize the probability of correct detection 

for a specified probability of false alarm. 

SOLUTION: Compute the likelihood ratio 

Likelihood of D(t) under Ej, 

“ Likelihood of b(t) under Hq 

Select Ej. when LR > threshold value. 

Select Hq when LR threshold value. 

As depicted in Fig. 3.7-1, the optimal decision rule can be 
efficiently implemented by processing the residual altimeter 
data D(t) with one matched filter and a threshold detector to 
test H.J, against Hq for all possible values of T. Once a de- 
tection is made (i.e., Ej. is selected), the maximum- likelihood 
estimates of the location T and the amplitude scale factor 
are easily computed from the matched-filter output. 
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Figure 3.7-1 Matched-Filter Detector 


Optimal Matched Filter - The optimal matched filter 
for long data sets is a convolution operator having the fre- 
quency response H(F) and the impulse response h(t) (e.g. , 
Ref. 18, pp. 325-329) 


h(t) 

H(F) 



gi2nFt 


M(-F) 


(3.7-1) 

(3.7-2) 


F = normalized frequency (cycles/sample) 

M(F) = Fourier transform of the ocean-current 
signature m(t) 

Sjjjj(F) = power spectrum of the residual altimetry N(t) 
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The Fourier transform of m(t) is defined as 


0» 

M(F) = m<t) 

tS-« 

•hi 

m(t> = I M(F) dF 


■L 


(3.7-3) 


(3.7-4) 


The power spectrum of N(t) is the Fourier transform of the 
autocorrelation sequence Rjjjj(t) 

Rjjj^(t) = E[N(t + t‘) N(f)l (3.7-5) 

oa 


c /p\ 


V' 


D / <- \ 

‘'NN''" 


^-i2nFt 


/O 7_A\ 
\ • r “ w / 


t=s-<» 


FIR Matched-Filter Algorithm - When the residual al- 
timeter noise model N(t) is autoregressive (AR) , the optimal 
matched filter can be implemented as a finite- impulse- response 
(FIR) filter. This means that the matched- filter impulse re- 
sponse h(t) has finite support, i.e., h(t) = 0 for t< and 

t > t . for finite t . and t . ' 

The AR model for the noise N(t) is a difference equa- 
tion of order p driven by zero-mean Gaussian white noise W(t) 

N(t) = N(t-l) + C 2 N(t-2) +...+ Cp N(t-p) + W(t) 

(3.7-7) 

t ** •*«— 1,0,1... 

Mean(W(t)) = 0; Variance (W(t)) = 
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The use of a Gaussian noise model is consistent with 
the observed Gaussian behavior of SEASAT-minus-geoid residual 
altimeter data described in Section 2.5. 

From linear-system theory (Refs. 17 and 18), the power 
spectrum of N(t) is 


" G(F)%(-F) 


(3.7-8) 


G(F) = 1 



i2nFk 


(3.7-9) 


From Eqs. 3.7-2 and 3.7-8, the matched filter frequency response 
may be expressed as 

H(F) = o"2g(F) G(-F) M(-F) (3.7-10) 


The Fourier transform of Eq. 3.7-10 yields the follow- 
ing expression for the impulse response h(t) of the matched 
filter (asterisks denote convolutions) 


h(t) 

= a“^g(t) * g(-t) * m(-t) 

(3.7-11) 

g(t) 

k 

= r G(F) dF 

“h 

(3.7-12) 

g(t) 

II 

o 

rt 

A 

O 

(3.7-13) 

g(0) 

= 1 

(3.7-14) 
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(3.7-15) 


g(p) = -Cp (3.7-16) 

g(t) = 0, t > p (3.7-17) 

Since the convolutions in Eq. 3.7-11 contain only a finite 
number of non-zero terms when m(t) has finite support, the 
impulse response h(t) also has finite support. Figures 3.7-2 
and 3.7-3 show examples of optimal impulse responses for cold- 
ring and Gulf-Stream signatures. 


Matched-Filter Performance Measures - The theoretical 
rms signal-to-noise ratio of the matched filter is defined as 


SNR 


Peak Filter Output Due to Signature m(t ) 
iRms Filter Output Due to Noise N(t) 


(3.7-18) 



Figure 3.7-2 Typical Impulse Response of Filter 

for Detecting Cold Rings 
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Figure 3.7-3 Typical Impulse Response for 

Gulf-Stream Detection 


The SNR of an optimal matched filter is an important 
parameter. As shovm in the Appendix, SNR may be computed with 
the foniiula 

SNR = .1 (3.7-19) 

1 j = -» 

where only a finite number of terms contribute. Vhe rms value 
of the noise in the filter output is numerically equal to SNR 
when the filter is optimized: 

Rms Noise in Filter Output = SNR (3.7-20) 
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The peak filter output value due to the signature m(t) is also 
expressible in terms of SNR: 

2 

Peak Filter Output Due to Signature mu ' = SNR (3.7-21) 

Equations 3.7-20 and 3.7-21 are derived in the Appendix. 

Since the SNR equals the rms value of the modeled 
noise in the filtered output, it is convenient to scale the 
filter output by the factor 1/SNR as shown in Fig. 3.7-1. 

This yields a convenient test statistic Y(t) that contains a 
random component havir a standard deviation of unity. 

The mean output frequency F^^, of the filter is a num- 
ber that measures the average rate at which the filter output 
changes sign 

F„ = Half the Average Rate of Zero Crossings /o 7 . 09 \ 
® of Noise in Filter Output ' 

As shown in the Appendix, the mean output frequency of the 
Gaussian noise in the filter output is given by 

Fm = ^ cos‘’^[X/SNR^] (cycles/sample) (3.7-23) 


where 


00 

X = V h(k) m(l-k) 
k=-w 


(3.7-24) 
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THRESHOLD DETECTION THEORY 


As depicted in Fig. 3.7-1 s the threshold detector 
compares the sequence of scaled test statistics Y(t) from the 
matched filter against a detection threshold. When Y(t) ex- 
ceeds the threshold, an alarm is said to occur. These alarms 
are classified into three categories; 

• Correct Detection caused by occurrences 
of modeled current signatures in the 
residual altimeter data 

• False Alarms caused by random excursions 
of the modeled noise in the residual 
altimeter data 

• Unmodeled Detections caused by unmodeled 
current signatures or unmodeled noise in 
the residual altimeter data. 

The statistics of correct detections and false alarms 
are computed for the specific ocean-current signature and the 
specific noise model for which the filter was designed. On 
the basis of these statistics, the expected average performance 
of the detector is predicted, and the detection threshold is 
adjusted for a desired tradeoff between detection probabilities 
and false-alarm rates. 

Formulas are listed below for computing the following 
performance statistics of the detector: the probability of 

false alarms; the average false-alarm rate, the maximum- 
likelihood estimates of signature location and amplitude (and 
their rras accuracies), and the probability of detecting a sig- 
nature with a prescribed amplitude. 

Probability of False Alarm - Let TH denote the detec- 
tion threshold and ?£ the probability of a false alarm. The 
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?£ is the probability that the noise component alone in Y(t) 
will exceed TH and is given by the standardized normal proba- 
bility distribution function: 

?£ = Prob{Y(t) > TH; noise alone} (3.8-1) 

00 

Pf = Q(TH) = f exp(-xV2) dx (3.8-2) 

^ Ih Tin 

The detection threshold that yields a prescribed probability 
of false alarm is 


TH = (3.8-3) 

Average False -Alarm Rate - The threshold detector 
processes data from individual tracks of altimeter data; it is 
often reasonable to set the detection threshold so that a speci- 
fied number of false alarms is expected to occur per unit dis- 
tance along the track (expressed in the units of alarms per 
data sample). This false-alarm rate (FAR) is computed as 

FAR = F_ exp(-TH^/2) (alarms/sample) 

“ (3.8-4) 

F = mean output frequency of modeled 
noise in filter output 

Equation 3.8-4 is derived from results in Ref. 17, p. 492. 

The expected number of false alarms (EN) along a track of data 
having N samples of the test statistic Y(t) is 


EN = N*FAR (3.8-5) 

Detection Threshold - It is recommended that the de- 
tection threshold TH be chosen to yield a specified false-alarm 
rate FAR in accordance with the equation 
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(3.8-6) 




TH = V-2 ln(FAR/Fjj^) 

Equation 3.8-6 is obtained by solving Eq. 3.8-4 for TH. 

Detected Signature Location and Amplitude - The 
maximum- likelihood estimate of the location of a detected 
signature is the value of t for which Y(t) achieves its local 
maximum value above the threshold TH. Let t^ denote this esti- 
mate of signature location. The Cramer-Rao (C-R) lower bound 
on the rms error of this estimate depends on the maximum scaled 
filter output Y(t^); 

C-R Lower Bound = YTt T ^ samples) (3.8-7) 

This lower bound on the variance of t^ can be re- 
expressed in units of kilometers 

«3 

C-R Lower Bound = 2t^F — y(t Y (3.8-8) 

*^m ■‘‘s '^o' 

V = Altimeter along- track 
velocity (m/s) 

f = Data sampling frequency 
(samples/sec) 

Nominal values for SEASAT altimeter data are v = 6700 m/s and 
f = 1.0 sample/sec. The false-alarm rate, can also be expres- 
sed in terms of distance: 


FAR (alarms/sample) 


V 10 


-6 


[FAR (alarms/megameter)] 


(3.8-9) 


The maximum-likelihood estimate of the signature am- 
plitude scale factor A^ is 
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(3.8-10) 


; - 

“ SNR 

and the rms (one-sigma) error in this estimate is 1/SNR, The 
maximum-likelihood estimate of the detected signature is then 

Probability of Detection - The probability of detecting I 

!l 

the model signature ra(t-t^) with a filter matched to m(t) and 
a detector operating with threshold TH is 

= Prob {Y(t^) > TH; m(t) + noise} (3.8-11) 

9 ii 

The peak filter output due to m(t-t^) is SNR (Eq. 3.7-21). | 

Therefore the peak scaled output is 

Y(t^) = SNR + (zero-mean unit-variance noise) (3.8-12) 

Equation 3.8-12 implies that Eq. 3.8-11 may be expressed in 
the equivalent form 

Pj = Prob = { noise > TH - SNR} (3.8-13) 

P^ = Q(TH - SNR) (3.8-14) ' 

For studying the ability of a matched filter to dis- 
criminate between different signatures, it is useful to compute 
the probability of detecting a signature m*(t-t^) when the i 

ji filter is matched to a different signature m(t). Let SNR be 

ij the rms signal- to-noise ratio of the filter with respect to 

I m(t) (Eq. 3.7-19). Let the peak output of the filter due to 

j! the other signature m'(t-t ) be denoted x'. If this peak out- 

I ( O o , 

p put occurs at time t^, then 

ii 

I 

i 

It 
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(3.8-15) 


7 


00 

X' = ^ h(k) in' (t' - k) 

Smmmm W 

1<=-C» 


The rms value of the noise in the 
<Eq. 3.7-20). Therefore, the rms 
(SNR') of the filter with respect 


filter output is SNR 
signal- to-noise ratio 
to signature m'(t-t^) is 


SNR' = x'/SNR 


(3.8-16) 


The probability of detecting signature ra'(t) with a detection 
threshold TH as in Fig. 3.7-1 is therefore 

Pd(m') = Q(TH - SNR') (3.8-17) 

In the special case where m'(t) = A m(t), the probability of 

s 

correct detection is 


Pd(As m) = Q(TH - A^-SNR) (3.8-18) 


3.9 SAMPLE PARAMETER VALUES 

Gulf-Stream Detection - Table 3.9-1 shows typical 
parameter values computed by the detection algorithm operating 
on residual SEASAT altimeter data from Rev. 234A. The table 
permits a comparison of parameter values for two Gulf-Stream 
model signatures that differ only in their heights. Since 
their widths are the same (100 km), both of these signatures 
can be optimally detected by a single matched filter. 
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TABLE 3.9-1 

SAMPLE PARAMETER VALUES FOR GULF-STREAM DETECTION 


SEASAT REV 23AA 

MATCHED-FILTER PARAMETERS 

THRESHOLD-DETECTOR PARAMETERS 

MODEL 

SIGNATURE 

SNR 

ni 

(cyc/Mm) 

C-R LOCATION 
ERROR BOUND 
(km) 

FAR 

(A/Mm) 

TH 

fd 

HEIGHT = 1.4 m 


■ 



10“^ 

5.5 

0.96 

WIDTH = 100 km 


4.1 

5 

lO"^ 

4.6 

1.00 



■ 



io“^ 

3.5 

1.00 

HEIGHT = 1.0 m 





1 

O 

m 

0.36 

WIDTH s 100 km 

5.2 

4.1 

f 

1 

..-4 

iU 


0.71 






2 

lO-"^ 

3.5 

0.95 


The larger 1.4-m signat-ure yields an rms signal-to- 
noise ratio (SNR) of 7.2 while for the smaller 1.0-m signature 
SNR = 5.2. The mean frequency (F„) of the noise in the filter 
outputs is 4.1 cycles/megameter for both filters, because 
is independent of the model signature amplitude. The Cramer- 
Rao (C-R) rms error bounds on locating the two signatures are 
5 km and 7 km; the smaller error corresponds to the larger- 
amplitude signature. 

The right side of Table 3.9-1 displays parameter val- 
ues for the threshold detector that compares the scaled output 
Y(t) of the matched filters with the threshold TH. Three 
false-alarm rates (FAR) are considered, ranging from 10 
alarm/megameter (A/Mm) to 10" A/Min. The detection thresholds 
(TH) range from 3.5 to 5.5 independently of the model signature 
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amplitude. The interpretation of the threshold values is sim 
pie because the test statistic Y(t) is scaled to have a vari- 
ance of unity. For example, with a threshold of 3.5, a false 
alarm occurs if the noise in Y(t) makes a positive excursion 
that exceeds 3.5 standard deviations from the mean (zero). 

The last column in the table lists the probabilities of de- 
tecting the signatures with the three false-alarm rates. The 
results show that the probability of detecting the smaller 
signature is sensitive to the false-alarm rate. The results 
in this table verify that typical Gulf-Stream signatures are 
expected to be reliably detected. 


Cold-F.ing Detection - Table 4.9-2 shows sample param- 
eter values computed by the detection algorithm operating on 
residual SEASAT altimeter data for Rev. 234A. This table com- 
pares parameter values for two cold=ring model signatures that 


differ in both width and amplitude. Because their widths are 
different, two different matched filters are required to detect 
these signatures optimally. 


The filter for the narrower (150-km) signature has 
the larger mean frequency (F^) and the lower Cramer-Rao loca- 
tion error bound. This occurs because the 150-km signature is 
more compact and leads to a higher signal-to-noise ratio than 
the 300-km signature. The 300-km signature is slightly more 
difficult to distinguish from the noise in the altimeter data 
as revealed by the values of SNR and the probabilities of de- 
tection (P^). 

A study of the parameter values for the threshold 
detector shows that the cold-ring signatures are detectable 
(P^ > 0.8) at false-alarm rates of FAR > 0.1 alarm/megameter 
(FAR > 1 alarm per 10,000 km of data). 
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TABLE 3.9-2 

SAMPLE PARAMETER VALUES FOR COLD- RING DETECTION 


SEASAT REV 234A 

1 

MATCHED-FILTER PARAMETERS 

THRESHOLD-DETECTOR PARAMETERS 

MODEL 

SIGNATURE 

SNR 

(cyc/Mra) 

C-R LOCATION 
ERROR BOUND 
(ko) 

FAR 

(A/Mm) 

TH 


DEPTH = 0.5 m 




0.01 

3.6 

0.61 

WIDTH = 150 km 

3.9 

6.7 

6 

0.1 

2.9 

0.84 





0.5 

2.3 

0.95 

DEPTH = 0.7 m 




0.01 

3.4 

0.49 

WIDTH = 300 km 

3.4 

3.0 

16 

0.1 

2.6 

0.77 





0.5 

1.9 

0.93 


3.10 ALGORITHM CAPABILITY FOR DISCRIMINATING BETWEEN 
TWO COLD RINGS 

Matched filters can be used to discriminate between 
cold-ring signatures having different widths. To study this 
discrimination capability, the probability of detection is 
plotted as a function of the false-alarm rate for signatures 
of different sizes. These plots are like the "operating char- 
acteristic" of the detector, except that the false-alarm rate 
is used in place of the false-alairm probability. 

As an example, consider a filter that is matched to a 
0.5-ra X 150-kra cold-ring signature (the signature is given by 
Eq. 3.3-1 with D 0.5 m and W = 150 km). For a typical track 
of residual SEASAT altimeter data in the western North Atlantic, 
the filter has an rms signal-to-noise ratio (SNR, Eq. 3.7-19) 
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of 3o9. The probability of detection (P^^, Eq. 3.8-14) for 
this SNR is plotted as function of the false-alarm rate (FAR, 
Eq. 3.8-4) in the upper part of Fig. 3.10-1, The different 
values of FAR correspond to different detection thresholds, 
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Figure 3.10-1 Filter Matched to 0.5-m x 150-km Signature 


When the 150-km filter is used to detect (suboptimally) 
a larger 0.7-m x 300-km cold-ring signature, the probabilities 
of detection are significantly reduced as indicated by the 
curve in the lower part of Fig. 3.10-1. (The SNR with the 
300-km signature is 2.2 as computed with Eq. 3.8-16, and the 
is computed with Eq. 3.8-17.) The vertical separation of 
the curves in Fig. 3.10-1 is a quantitative measure of the 
discrimination capability of the 150-km filter as a function 
of FAR. 


Figure 3.10-2 is like Fig. 3.10-1, except the matched 
filter is optimized here to detect the wider 300-km signature. 
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Because the filter is optimized for the wider signature, the 
SNR (3.4) is larger than it is in Fig. 3.10-1, where SNR = 2.2 * 
The performance of the 300-km filter for detecting (suboptiraally) 
the 150-km signature is depicted by the curve in the lower 
part of Fig. 3.10-2. A comparison of Figs, 3.10-1 and 3.10-2 
shows that the wider 300-km filter is less discriminating and 
more robust than the narrower ISO-km filter, as quantified by • 
the smaller separation of the curves in Fig. 3.10-2, 
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Figure 3.10-2 Filter Matched to 0.7-m x 300“km Signature 


The results presented here for detecting cold rings 
are considered typical and may be summarized as follows: 


• The filter matched to the narrower sig- 
nature is the more discriminating 

• The filter matched to the wider signa- 
ture is the more robust. 
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3.11 SUMMARY OF OCEAR-CURRENT DETECTION STUDY 


An ocean-current detection algorithm for processing 
single tracks of residual satellite altimeter data has been 
developed and verified. The algorithm is data-adaptive and is 
based on autoregressive statistical modeling, matched detection 
filters, and models for generic altimetric signatures caused 
by boundary and ring currents. 

The algorithm performance was verified for Gulf-Stream 
and cold-ring detections using SEASAT altimeter data, the Marsh- 
Chang geoid, and surface truth data in the western North Atlantic. 
The algorithm performance was quantified in terms of 

• False-Alarro Rates 

• Probabilities of Detection 

• Location Estimation Errors 

• Signal-to-Noise Ratios 

• Mean Output Frequencies. 

The study results indicate that the automatic detec- 
tion of Gulf-Stream and cold-ring signatures is feasible with 
a processing time (IBM Model A341 computer) of approximately 
0.5 sec per data track using one matched filter. The algorithm 
is also applicable to warm-ring detection. 

The results of verifying the algorithm with SEASAT 
data lead to the following conclusions: 

• The algorithm performance is insensitive 
to localized data errors 

• Filters matched to narrower signatures 
discriminate against wider signatures 
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• Filters matched to wider signatures are 
robust; i.e., they are especially insen- 
sitive to localized data errors, and 
they are able to detect (suboptimally) 
narrower signatures of similar shape. 


Formal specifications for the detection algorithm are 
given in Ref. 19. 


3.12 AREAS FOR FUTURE STUDY IN OCEAN-CURRENT DETECTION 


The following three areas are suggested for future 

study. 


• Quantify the sensitivity of the algorithm 
^rformance to geoid estimation* error . 
Although the algorithm is data adaptive 
and can therefpre be used with any geoid 
estimates, the sensitivity of the algo- 
rithm performance to the accuracy of the 
geoid estimates remains to be studied. 

• Develop a recursive implementation for 
the matched filters . Because the altim- 
eter data are modeled autoregressively , 
the matched filters can be rigorously 
implemented with a recursive algorithm, 
as opposed to the straight convolutions 
used in the present algorithm. A recur- 
sive implementation will reduce the time 
to filter the data by a factor of five 
and will be useful for applications to 
multi-signature detectors for ring currents. 

• Process all relevant SEASAT data in the 
western North Atlantic . This will pro- 
vide additional statistics describing 
the algorithm performance in an opera- 
tional setting. The western North 
Atlantic has the distinction of being an 
oceanographically interesting region for 
which geoid estimates and surface truth 
data are available. 
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A. GEQSTROPHIC-VELOCITY ESTIMATION ALGORITHM 

A . 1 INTRODUCTION 

A nearly geostrophic boundary current, such as the 
Gulf Stream in the western North Atlantic, produces a charac- 
teristic signature in satellite altimeter data where the sat- 
ellite subtrack intersects the current. The along-track slope 
of this signature is proportional to the cross- track component 
of the geostrophic current velocity as ''ixpressed by the geo- 
strophic relation. Therefore, the cross- track geostrophic 
velocity can be inferred from estimates of the signature slope. 
This chapter describes a velocity estimation algorithm based 
on estimating signature slopes via matched filtering. 

The chapter is organized as follows. Section A. 2 
describes an automated data-adaptive algorithm for computing 
estimates of geostrophic velocity. The algorithm uses the 
matched-filtering detection algorithm described in Chapter 3. 
Verifications of the velocity estimation algorithm with SEASAT 
altimeter data are discussed in Section A. 3. A summary and 
conclusions are presented in Section A. A, and suggestions for 
further study are offered in Section A. 5. 


A . 2 ALGORITHM DESCRIPTION 

The velocity estimation algorithm processes single 
tracks of residual satellite altimeter data (sea surface height 
minus geoid estimates). Each track of data is assumed to have 
been previously processed with a single matched filter that 
detected a boundary-current signature in the data. 
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The outputs of the algorithm are the estimated along- 
track location of the maximum geostrophic velocity, the esti- 
mated maximum cross- track geostrophic velocity, and rms 
accuracy measures for the location and velocity estimates. 

The algorithm uses several matched filters to esti- 
mate the parameters and A_ in the following hyperbolic- 

s s 

tangent model signature (which is discussed in Section 3,3): 

H(x) = - (A/2) tanh(3 x/W„) (4.2-1) 

• s s 

H(x) = Sea surface height 

X = Position along subtrack 

A = Amplitude of model signature (arbitrary 

initial value used to design matched filter) 

A = Amplitude scale factor (chosen for best fit 
with data) 

W = Signature width (90 percent height change). 

The cross-track geostrophic velocity is estimated 
from the geostrophic relation applied to the slope of the model 
signature. By using equations in Section 3.3, the cross- track 
velocity implied by the model signature can be written as 

= ^s r f -| - sech2(3 x/W ) (4.2-2) 

s 

V^(x) = Cross- track component of geostrophic velocity 
X = Position along subtrack 
f = 2fisin(j) = Coriolis parameter 
Q = Earth's rotational velocity 
(j> = North latitude 
g = Acceleration of gravity. 
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The maximum cross-track geostrophic velocity is 

Nmax = V,(0) = (4.2-3) 

The uncertainty in the velocity estimate is estimated from the 

rms uncertainty in the estimate of parameter A_. As discussed 

s 

in the following, the algorithm consists of five steps. 

Step 1 - The residual altimeter data are processed 
with five matched filters, each optimized for detecting a 
different-width signature. For example, signature widths of 
50, 60, 75, 100 and 150 km are reasonable for the Gulf Stream 
in the western North Atlantic. 

Step 2 - The most likely width of the signature in 
the altimeter data is estimated by comparing the scaled out- 
puts (Y(t) in Fig. 3.7-1)) from each of the matched filters to 
determine that filter which yields the largest response to the 
signature. This filter is said to be best-matched to the sig- 
nature in the data. The signature width of this filter is 

A 

selected as the best estimate for the width of the signature 
in the data. 

A 

Step 3 - The maximum- likelihood estimate A^ of the 
a.:^)lirude scale factor is computed from the observed maximum 
output Y(t ) of the best-matched filter in accordance with the 
theory discussed in Section 3.8: 

K - Y(t )/SNR (4.2-4) 

SNR = Theoretical rms signal- to-noise ratio of the 
best-matched filter 
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SNR is computed by using Eq. 3.7-19. Since Y(t^) is modeled 
as a sample from a random process with unit variance, the rms 

A A 

uncertainty 6A_ in the estimate A_ is given as 
s s 

6A_ = 1/SNR = Standard deviation of A_ (4.2-5) 

s s 

The estimated maximum cross-track geostrophic velocity 

A A 

V^(0) and its rms uncertainty 6V^(0) are given as 

V (0) = A_ (4.2-6) 

® 2 f Wg 


6V^(0) = 6^ ^-^-4- (4.2-7) 

^ ® 2 f Wg 

Equation 4.2-7 is a useful working estimate of veloc- 
ity estimation accuracy. It is recognized that a more compli- 
cated error analysis is possible, which accounts for the error 

A A A A 

6W^ in and the correlation between 6W„ and 6A„. Since no 
s s s s . 

stochastic model is available to support such an analysis, the 

0 

simpler estimate in Eq. 4.2-7 is recommended. 

Step 4 - The maximum- likelihood estimate of the loca- 
tion of the maximum geostrophic velocity is that time t^ at 
which the best-matched filter achieves its maximum scaled out- 
put Y(t^). The rms uncertainty 6t^ of this estimate is given 
by the Cramer-Rao lower bound (Eq. 3.8-7): 

6to = in Y(t^) (4.2-8) 

where the mean frequency is computed by using Eq. 3.7-23. 



t; 


4.3 DEMONSTRATION OF GULF-STREAM VELOCITY ESTIMATION 

The performance of the velocity estimation algorithm 
is verified by using it to estimate Gulf-Stream velocities. 

Four examples are presented, which use the same tracks of SEA- 
SAT residual altimeter data as the demonstrations of the de- 
tection algorithm in Sections 3.4 and 3.5. The time required 
to process a track of data containing one boundary-current 
signature is 1.5 seconds on an IBM Model 4341 computer. (Only 
the data in the vicinity of the detected boundary current are 
used for geostrophic velocity estimation.) 

Example 1 - Figure 4.3-1 shows two plots of sea sur- 
face height as a function of position along the SEASAT sub- 
track for Rev. 234A. The solid curve is the residual altimeter 
data that were used as input to the velocity estimation algo- 
rithm described in Section 4.2. The dashed curve is the tanh 

I 

model signature that the algorithm fitted to these data. The 
vertical positioning of the model signature has no effect on 
the estimate of geostrophic velocity, which depends only on 
the slope of the signature. Therefore, the algorithm does not 
estimate the vertical positioning of the signature; the verti- 
cal position of the tanh signature is manually selected in 
these examples to provide easy visual comparisons with the 
SEASAT -data. In contrast, the vertical and horizontal scales 
and the horizontal position of the tanh signature are each 
important for estimating the current velocity; these features 
of the signature are determined automatically by the algorithm. 
The fit of the tanh signature to the SEASAT data looks reason- 
able in the Gulf-Stream transition region centered on data 
sample number 233. The maximum cross- track component of the 
geostrophic velocity is 1.54 m/s (Eq. 4.2-6) with an estimated 
rms accuracy (Eq. 4.2-7) of 0.22 m/s (14 percent of the velocity 
estimate) . 
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Figure 4.3-1 Gulf-Stream Signatures, Rev. 234A 


Example 2 - Figure 4.3-2 shows data for SEASAT Rev. 
277A. The estimation algorithm provides a reasonable fit to 
the altimeter data in the transition region centered on sample 
number 233, where the algorithm estimates that the cross-track 
velocity achieves a maximum value of 2.82 m/s. The rms accu- 
racy is estimated by the algorithm to be 0.70 m/s (25 percent 
of the velocity estimate). This estimate of rms uncertainty 
may be too large, because in this example the algorithm used 
data (for noise modeling) that contained significant data errors 
due to the island of Bermuda as depicted in Figure, 3.4-4. 

The deviation between the tanh signature and the altim- 
eter data on the right of Fig. 4.3-2 is caused by warm ring U, 
whose existence is verified by independent oceanographic data 
(Ref. 6). Fortunately, the estimation of maximum geostrophic 
velocities requires only that the sea surface slope be estimated 
in the relatively narrow transition zone where the geostrophic 
flow is located. 
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Figure 4.3-2 Gulf-Stream Signatures, Rev. 277A 


Examples 3 and 4 - The model signatures selected by 
the estimation algorithm for SEASAT Revs. 478A and 572D are 
shown in Figs. 4.3-3 and 4.3-4. The tanh signatures provide 
reasonable fits with the altimeter data in the transition 
regions of the Gulf Stream. The results in Fig. 4.3-4 demon- 
strate that the fit of the model signatures to the maximum- 
slope of the altimeter data is insensitive to nearby ring cur- 
rents. The estimated maximum cross- track geostrophic velocity 
for Rev. 478A is 1.78 m/s with an estimated rms accuracy of 
0.36 m/s (20 percent of the velocity estimate). For Rev. 572D 
the estimated maximum velocity is 1.96 m/s with an estimted 
xrms accuracy of 0.26 m/s (13 percent ol the velocity estimate). 
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Figure 4.3-4 Gulf-Stream Signatures, Rev. 572D 
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A. 4 SUMMARY AND CONCLUSIONS OF GEOSTROPHIC-VELOCITY 
ESTIMATION STUDY 

Sununary - An algorithm has been developed to estimate 
the velocities of geostrophic surface currents by the analysis 
of single tracks of satellite altimeter data. The algorithm 
performance is verified by four examples of Gulf-Stream veloc- 
ity estimation using SEASAT altimeter data and the Marsh-Chang 
geoid. 


The algorithm utilizes the data-adaptive matched- filter 
algorithm developed for ocean-current detection and automatically 
performs the following analysis of the altimeter data. A bank 
of five matched filters estimates the location and the best 
width and amplitude for a model signature of the sea-surface 
height where the satellite subtrack crosses the geostrophic 
current. The cross- track component of the geostrophic current 
velocity is estimated from the slope of the model signature by 
using the geostrophic relation. The rms accuracy of the veloc- 
ity estimate is computed from the theoretical rms uncertainty 
of a maximum- likelihood estimate of the signature amplitude. 


Conclusions - The results of this study lead to the 
following conclusions; 


• The automated estimation of cross-track 
geostrophic velocities from satellite 
altimeter data is feasible for boundary 
currents such as the Gulf Stream 

• The velocity estimation algorithm is a 
simple extension of the data-adaptive 
matched-filter algorithm for detecting 
and locating boundary and ring currents 

• Processing time is about 1.5 sec per 
data track using an IBM Model 4341 com- 
puter; 0.5 sec for detecting the current 
signature with one matched filter, and 
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1 sec for estimating the geostrophic 
current velocity with four additional 
filters operating on the data in the 
vicinity of the current signature 

• The performance of the velocity estima- 
tion algorithm is insensitive to isolated 
data errors and to the occurrence of 
warm-ring and cold-ring signatures in 
the altimeter data. 


Formal specifications for the velocity estimation 
algorithm are given in Ref. 19. 


A. 5 SUBJECTS FOR FURTHER STUDY 


Based on the results of this investigation, the fol- 
lowing subjects are recommended for further study. 


« Determine the best choice of model sig - 
natures 7 The use of five model signa- 
tures in this study (50, 60, 75, 100, 
and 150 km) yields reasonable results, 
but other choices may be preferable in 
terms of the trade-offs between accuracy 
and computational complexity. 

• Quantify the sensitivity of the velocity 
e stimates to the accuracy of the geo 15 
estimates . ATthough the algorithm is 
data-adaptive and can be used with any 
geoid estimates, the sensitivity of the 
velocity estimates to geoid error has 
not been studied. 

• Compare current velocities estimated by 
the algorithm with surface tHith data . 
Because no surface truth data on Gulf- 
Stream velocities were available for the 
subtracks used in this study, a direct 
comparison of satellite-derived velocities 
with surface data has not been possible. 
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• Develop a faster recursive algorithm for 
Implementing the matched filters in the 
algorithm . Sucfi recursive algorithms 
are feasible , because the altimeter data 
noise model is autoregressive. A recur- 
sive algorithm would permit larger numbers 
of model signatures to be used for in- 
creased accuracy. Alternatively, recur- 
sive filters could be used to reduce the 
computational requirements. 
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SUMMARY. CONCLUSIONS, AND RECOMMENDATIONS 


5 . 1 SUMMARY 

This report documents the study of the oceanographic 
information content of SF.ASAT altimeter data and describes the 
development of NOSS algorithms for ocean-current mapping and 
their verification with SEASAT altimeter data. The inputs to 
the algorithms are individual tracks of residual satellite 
radar altimeter data from which estimated geoid profiles have 
been subtracted. The algorithms are based on the fact that 
cold-core and warm-core current rings and boundary currents 
can be detected by identifying the occurrence of characteristic 
sea-surface height signatures in the residual altimeter data. 

In the case of nearly geostrophic boundary currents, the cross- 
track component of the current velocity can be inferred by 
estimating the along-track sea-surface slope, from the altimeter 
data and then using the geostrophic equation to compute the 
velocity. 


In the study of the oceanographic information content 
of SEASAT altimeter d?ta, autoregressive modeling techniques 
were used to analyze the statistical behavior of residual al- 
timeter data. The analysis results show that residual altim- 
eter data have the statistical properties of Gaussian random 
noise with a correlation structure that varies from track to 
track. The results of the study were used as a basis for the 
development of the NOSS algorithms. 

The NOSS algorithms are based on optimal matched fil- 
ters, which are used to detect, locate, and estimate the ampli- 
tudes of generic current signatures in the residual altimeter 
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data. The algorithms automatically analyze each track of resid- 
ual altimeter data and compute an optimal autoregressive model 
for the noise signal in the data. Using this noise model, 
together with a parametric model for the deterministic ocean- 
current signature that is to be detected, the algorithm designs 
a statistically optimal matched-filter detector for discriminat- 
ing between the noise and the signature. The detector is op- 
timal in the sense that the probability of detecting the ocean- 
current signature is maximized for a specified probability of 
false alarm (a false alarm occurs when the random noise excur- 
sions in the altimeter data masquerade as a current signature 
and cause a false detection). The algorithm adjusts the sensi- 
tivity of the detector to achieve a specified average false- 
alarm rate (e.g., 1 false alarm per 10,000 km of data along 
the satellite subtrack). The algorithm perfomance was demon^^ 
strated with SEASAT altimeter data and Marsh-Chang geoid esti- 
mates in the western North Atlantic, where known boundary- 
current and ring-current signatures occur in the altimeter 
data. 


The algorithm for estimating the geostrophic veloci- 
ties of boundary currents employs a bank of five matched- filter 
detectors; each filter is matched to a different width for the 
current signature. The algorithm determines that signature 
width which is most probable (given the available altimeter 
data) and computes a maximum- likelihood estimate of the current 
signature amplitude. From this information, the algorithm 
estimates the maximum along-track slope of the sea surface and 
uses the geostrophic equation to compute the estimated cross- 
track component of geostrophic velocity. The ms accuracy of 
the velocity estimate is also computed by using the Cramer-Rao 
lower bound on the variance of the estimated signature ampli- 
tude. The algorithm perfomance was demonstrated with SEASAT 
data containing Gulf-Stream signatures. 


5-2 







The formal specifications for the NOSS algorithms are 
documented in Ref. 19. 


5.2 CONCLUSIONS 


The primary conclusions of this study are summarized 
in the following. More detailed discussions are provided in 
Section 2.8, 3.6, 3.11, and 4.4. 


• The residual SEASAT altimeter data (SSH 
minus geoid) analyzed in this study has 
statistical properties of correlated 
Gaussian random noise. No systematic 
departure from Gaussian behavior was 
observed. 

• By using a unique Kalman smoother to 
process each track of data, the time- 
varying mesoscale SSH waveforms in sets 
of difference data from nearly repeating 
SEASAT tracks can be estimated with an 
average theoretical rms accuracy of 2 cm. 
The results of using such smoothers on 
data from the region south of Iceland 
indicate that statistically significant 
mesoscale variations (larger than 6 cm) 
occur during several of the 3-day inter- 
vals between Revs. 

• Verification tests with residual SEASAT 
altimeter data from the western North 
Atlantic indicate that the Gulf Stream 
is reliably detected and located with 
the NOSS detection algorithm developed 
in this study. The required computing 
time is less than 0.5 sec per track using 
an IBM 4341 computer. 

• The feasibility of using the NOSS detec- 
tion algorithm to detect and locate cold- 
core current rings in the western North 
Atlantic was verified with residual SEASAT 
altimeter data. To increase the probabil- 
ity of detection and to discriminate be- 
tween rings of different widths, it is 
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recommended that each track of data be 
processed by the algorithm with two or 
three matched filters optimized for 
different-width current signatures. 

♦ Verification tests with SEASAT data indi- 
cate that the automated estimation of 
cross-track geostrophic velocities from 
residual satellite altimetry is feasible 
for boundary currents such as the Gulf 
Stream. The required processing time 
with Che NOSS algorithm is about 1.5 sec 
per data track with ar. IBM 4341 computer. 

• Both the current-detection and velocity- 
estimation algorithms yield outputs that 
are insensitive to isolated data errors. 


5.3 SUBJECTS FOR FURTHER STUDY 


Based on the results of this investigation, the 
following subjects are recommended for further study. More 
detailed suggestions are provided in Sections 2.8, 3.12, and 
4.5. 


e Develop statistical models for the tem- 
poral variations in SEASAT significant- 
wave-height data and their correlations 
with the observed mesoscale variability 
over 3-day intervals in SEASAT repeat- 
track SSH data. 

• For both current detection and velocity 
estimation, quantify the sensitivity of 
the NOSS algorithm performance to geoid 
estimation error. 

• Develop a recursive implementation for 
the matched filter in the NOSS algorithm 
to provide increased computational 
efficiency . 

• Process all relevant SEASAT altimeter 
data in the western North Atlantic to 
provide additional statistics describing 
NOSS algorithm performance. 
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APPENDIX A 

DERIVATIONS OF MATCHED- FILTER EQUATIONS 


This appendix brings together in a consistent notation 
the derivations of .several key results in matched- filter theory. 
Previously published derivations span a variety of mathematical 
notations, are scattered in the literature, and do not address 
the discrete- time theory in detail. 


Derivation of Eg. 3.7-20 - The frequency response of 
the optimal matched filter is (e.g.. Ref. 18, pp. 325-329) 


H(F) 


- M(-F) 
■ 
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= complex conjugate) 


Snn(F) = Power spectrum of noise in data. 


Let x(t) denote the output of the filter at time t. 
When the data consist of noise only, the mean-square value of 
X(t) is 
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Since M(-F) = M*(.F) = complex conjugate of M(F), 


E(x^(t)l 

E(X^(t)I 

E[x^(t)l 


I 


£ 

J 


M(F) dF 


H(F) M(F) dF 


H(F) M*(-F) dF 


00 


h(k) m(-k) 


k=-oo 


where 


(A-4) 

(A-5) 

(A-6) 

(A-7) 


h(k) 

m(k) 



^i2nFk rxr 
e DF 


gi2nFk 


(A-8) 

<A-9) 


Therefore, the rms value of the noise in the output of the 
optimal matched filter is 


RMS 
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Derivation of Eqs. 3.7-20 and 3.7-21 - The maximum 
output from the filter specified by Eq. A-1 occurs at 

lUo^ 

time t=0 when the input is ra(t) 
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where 
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But X(F) = H(F) M(F), therefore 
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Comparing Eq. A-13 with Eq. A-5 shows that 


X = Efv^ftn 
*max 


(A-14) 


This result can be rewritten as follows by using Eq. A-10 


X = RMS' 
^max 


(A-15) 


The rms signal- to-noise ratio (SNR) is defined as 


SNR = 


“max 




(A-16) 


Substituting Eqs. A-15 and A-10 in Eq. A-16 yields 


SNR = ^ = RMS 


(A-17) 


which verifies Eq. 3.7-20. Equation 3.7-21 is verified by 
using Eq. A-17 to rewrite Eq. A-15 as 


A- 3 





Oi; - ^ ■' " ‘ 

X„,^v = SNR^ (A-18) 

max 

Derivation of Eg. 3.7-23 - Let U(t) denote the unit 
step function defined as 


U(t) = 1 for t > 0 
U(t) = 0 for t < 0 


(A-19) 


The average number C(N) of zero crossings in N samples of a 
stationary time series {X^} is 


C(N) 



N 


k=l 


U(-X 


k ^k+l' 


(A-20) 


The mean frequency 


expected value of C(N) 


V nf 

‘m 


{Xj^} is defined as half the 


F„ = I EIC(N)) 


_ 1 
" 2n 


N 


k=l 



(A-21) 


Because {Xj^} is stationary, the expectation in Eq. A-21 is inde- 
pendent of k. Therefore, Eq. A-21 may be written as 

Fm = I (A-22) 

But E[D(-Xj^ ” Probability that X^ ^k+1 

expressed as 

= I Prob Xj^^j < 0} (A-23) 

Let 


r 


^xx<l>/^xx 


( 0 ) 


A-A 


(A-24) 





OF Pr-..J o-^-T7 

where is the autocorrelation sequence of 


® iXj+k Xji 

The analysis in Ref. 17, p. 199, shows that 

Prob(Xj^ Xj^^^ < 0) = B/n 

B = cos"^(r) 

From Eqs. A-23, A-26, and A-27, it follows that 

F - cos“^(r) 
m “ 2n 


(A-25) 


(A-26) 

(A-27) 


(A-28) 


The last step in this derivation is to establish a 
formula for parameter r, which is valid when {Xj^} is the out- 
put of an optimal matched filter with frequency response H(F) 
(Eq. A-1) that is driven by zero-mean Gaussian noise with power 
spectrum Sjjj^(F). The power spectrum S^(F) of {Xj^} is 

= |H(F)|2 Sjjjj(F) (A-29) 

In view of Eq. A-1, Eq. A-29 may be written as 




IM(-F)|^ 


(A-30) 


But |M(-F)|^ = M(-F) M(F) because m(k) is real. Therefore 

= ipr 

Sj^,^(F) = H(F) M(F) 

A- 5 


(A-32) 





OF POOK 

Taking the inverse Fourier transform of Eq. A-32 yields the 
following expression for the autocorrelation sequence of the 
output noise 




E 

j=-» 


h(j) m(k-j) 


(A-33) 


h(k) = Impulse response of filter 
m(k) = Signature to be detected 


With this result, it follows from Eqs. A-10 and A-17 that 

= SNr2 <A-34) 

Let X = R„„(l), then 

XX 

OB 

j=-“ 


Combining Eqs. A-34, A-35, A-24, and A-28 yields the desired 
result: 

^m ■ ^ cos’^(x/SNR^) (A-36) 
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